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Abstract 

Coupled distinct arrays of nonlinear oscillators have been shown to have a regime of high fre- 
quency, or ultra-harmonic, oscillations that are at multiples of the natural frequency of individual 
oscillators. The coupled array architectures generate an in-phase high-frequency state by coupling 
with an array in an anti-phase state. The underlying mechanism for the creation and stability of 
the ultra-harmonic oscillations is analyzed. A class of inter-array coupling is shown to create a 
stable, in-phase oscillation having frequency that increases linearly with the number of oscillators, 
but with an amplitude that stays fairly constant. The analysis of the theory is illustrated by 
numerical simulation of coupled arrays of Stuart-Landau limit cycle oscillators. 

PACS numbers: 05.45.-a, 05.45.Xt 
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I. INTRODUCTION 



The dynamics occurring in coupled oscillators has been of much interest in both basic 
and applied science. Since the original early observations of Huygen's of two synchronized 
coupled nonlinear clocks (See the appendix for a nice description in theories of coherent 



motion in 
numerous 
single statistical 




imit cycle and phase oscillator arrays having different couplings have become 
In most cases studied, however, the arrays are assumed to have a 
coup llng such as.oba! coupBng |f Q, nearest neighbor and 

ring coupling j9|,jl0], and long range coupling to list just a few. 

In many of the examples of coherent motion studied, the synchronized in-phase state and 
incoherent splay or anti-phase state have both been observed. In coupled semi-conductor 
lasers, it is known that nature prefers the incoherent states , while in Josephson junction 
arrays which are globally coupled, both in-phase and incoherent states may co-exist 13|, [14 1 . 
The stability of the in-phase state in most applications has been important in order to 
maintain frequency control and/or power output in the arrays, such as in coupled lasers and 
Josephson's junctions. 

Recently, interacting arrays of limit cycle oscillators possessing different architectures have 

n 

been proposed to stabilize a high frequency in-phase state [15] . While numerical simulations 
exhibit stable in-phase states over a particular range of parameters, no analytical work has 
been done on the system, and the mechanism behind the onset of stable ultraharmonic 
oscillations has not been explained. The goal of this paper is to explain the mechanism 
behind the onset of stable ultraharmonics, obtain analytic results for the critical value of 
the coupling constant required for this onset, as well as to explore an alternative coupling 
scheme previously not considered. The paper is organized as follows: In Section II, the 
equations of motion are introduced and solved exactly, in the absence of global coupling 
between the two arrays. In Section III, the bifurcation value of the global coupling constant 
above which stable ultra-harmonic oscillations set in is calculated and shown to depend of 
frequency, amplitude and total number of oscillators. The physical mechanism behind the 
onset of stable ultra-harmonics is explained, as resulting from a large constant component in 
the amplitude-dependent global coupling, that determines the bifurcation value. Above the 
bifurcation value, the dynamics of one of the arrays become a function of the other, resulting 
in ultraharmonics. Averaging theory is applied in Section IV to explain the equivalence of 
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the bifurcation diagrams between the averaged and the full system and to estimate the 
amplitude of ultraharmonic oscillations. Section V derives an alternative type of coupling 
that while more complicated in form, has the advantage of creating a higher amplitude 
ultra-harmonic oscillations, with an amplitude that stays fairly constant as the number of 
oscillators in the array increases. Section VI concludes and summarizes our results. 

II. BASIC EQUATIONS AND DYNAMICS 

The internal dynamics of the uncoupled Stuart-Landau oscillators in normal form is: 

z = (a + iuS)z — \z\ 2 z (1) 

where z is a complex variable. Equation (|1J) can be easily solved by transforming into polar 
variables, and results in a limit cycle oscillator with frequency equal to uj and steady-state 
amplitude of oscillation equal to yfa. The system considered here consists of two arrays 
of coupled Stuart-Landau oscillators, {Xj} and {Yj}, where the index, j stands for each 
individual oscillator. Throughout the rest of the paper, the {Xj} and {Yj} will be referred 
to simply as X and Y arrays. In the absence of coupling, the equations of motion of each 
oscillator within the arrays is given by Eq. We start by considering two coupled arrays, 
as initially introduced in Q|. The difference between the two arrays is that X oscillators 
possess anti-phase diffusive coupling between oscillators within the array, while the Y array 
has in-phase diffusive coupling within the array. The two arrays are globally coupled to each 
other. The whole system is described by the following set of equations: 

N 

Xj = (a x + iuo)Xj - \Xj\ 2 Xj + c x (X,_! + X j+1 - 2Xj) + c yx ^ \Y k \ (2) 

k=l 

N 

Yj = {a y + iu)Yj - \Y 3 \% + cy + Y j+1 - 2Yj) + c xy ^ \X k \ (3) 

k=l 

A schematic diagram of the architecture of the coupled arrays in Eqs. (J2j and is shown 
in Figure In the case of the X array oscillators, the diffusive coupling constant of 
c x < between oscillators within the array leads to an anti-phase solution. For the F-array 
oscillators, where c y > 0, the in-phase state is the only stable solution. Thus in steady-state, 
all oscillators in the F-array tend to have the same dynamics, so that the diffusive coupling 
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term in Eq. Q drops out. The V-array oscillators will all oscillate in-phase in steady-state, 
which is a direct outcome of the diffusive coupling term, c y > 0, so that the in-phase state 
is the only stable state. From Eqs. (j2J) and (jHJ), it is clear that global coupling depends 
only on the amplitude of oscillators in the two arrays, \Y^\ and \X^\. If this amplitude is 
approximately conserved, then the main contribution from the global coupling term will be 
a constant component. In Section V, an alternative coupling scheme that results in a higher 
oscillatory component in the global coupling, leading to higher intensity ultraharmonics, will 
be considered. 

Dynamics in the absence of inter-array coupling 

In the absence of inter-array couping, c xy = c yx = 0, an exact steady-state solution can be 
obtained for the two arrays. Since in steady-state, all the F-array oscillators move in-phase, 
the diffusive coupling term drops out of the Eq. (JHJ), and the solution reduces to that of a 
single limit cycle oscillator, so that the index j on the Y array will be dropped throughout 
the rest of the paper. Unlike the Y array oscillators, the diffusive coupling term in the 
X array is negative, so that the oscillators tend to be maximally out-of-phase with their 
nearest neighbors. For a ring of diffusively coupled oscillators, this results in two solutions, 
depending on whether N is even or odd: 

ACTi = * (4) 

A0°f +1 = 7T + 7T/N (5) 

where A^I^S and A(f)j^ +1 is the phase difference between the jth and the j + 1th oscillator 
in the X-array for TV-even and iV-odd, respectively. Since all the oscillators in the X-array 
are identical and differ only by a phase-shift, we now have an expression for the steady-state 
dynamics of j + 1 and j — 1 oscillators as a function of the jth oscillator: 

Xj-x = e-^'^Xj X j+1 = e iA ^ +1 Xj (6) 

Substituting Eq. (0) into Eq. (J2J), using Eqs. (jlj) and (0), and setting c yx = 0, we can 
now solve for the amplitude of the X-array oscillators in the absence of global inter-array 
coupling: 
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for even N, and 

A^° dd = [a x - 2c x (1 + cos^/N))} 1 ' 2 (8) 

for odd N. The frequency of X array oscillators is the same as in the Y array, the limit 
cycle frequency, to. For a x = a y however, since c x < 0, the amplitude of oscillation is higher 

1/2 

than that of the in-phase y-array, which oscillates at amplitude a y . As will be shown in 
the following section, this difference in amplitude is important to the creation of the stable 
ultraharmonics and affects the range of the global coupling constant necessary for stable 
ultraharmonics, in the case of symmetric coupling. 

From Eqs. (j3J) and (JSJ), it is clear that for iV-even, the X array splits into two identical 
subgroups containing N/2 oscillators, where the two subgroups are out-of-phase with each 
other by 7T. For iV-odd, the X-array forms a traveling wave, with phases distributed over 
the unit circle in increments of 2ir/N. This can be seen by using Eq. (JSJ), which gives 2ir/N 
as the smallest phase difference (between jth and j + 2nd oscillators). Since for generation 
of ultraharmonics at frequency Nu, the X array must form a traveling wave, with phase 
increments of 27r/N, the rest of the paper will focus on iV-odd arrays. The reason why this 
2n/N distribution of phases is needed for the creation of ultraharmonics at frequency Nuo 
will be explained in the following section. Other types of inter-array coupling schemes could 
be considered that would result in a different phase distribution. 



III. THE EFFECT OF GLOBAL COUPLING AND CREATION OF STABLE UL- 
TRAHARMONICS VIA A BIFURCATION 

It was shown in the previous section that in the absence of global coupling, c yx = c xy = 0, 
the steady state dynamics of the two arrays are that of a simple harmonic oscillator, with 

1 /2 

amplitudes given by a y and Eq. (jSJ), for the Y and X arrays, respectively. Since the global 
coupling term is a function of amplitude only (see Eqs. (J2J and (JHJ)), we would expect the 
global coupling to introduce a constant component into Eqs. (J2J) and Q. In fact, if the 
coupling is one-way, that is if we set c yx = 0, then the global coupling term only contributes 



a constant component and no u 
previously made observation 



traharmonics can occur in the Y array. This explains the 
that mutual coupling between the two arrays is required 
to induce ultraharmonic oscillations. It follows that the global coupling term in the X array 
must create a perturbation in the amplitude of the X oscillators, so that it is no longer 
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a conserved quantity. To understand how this happens, we examine the weakly coupled 
case, c yx ,c xy < 1/ (oix(yNj. Then, in the absence of resonant interactions, the steady-state 
dynamics of the two arrays is given to lowest order by: 

Xj = (a x + iuj)Xj - \Xj\ 2 Xj + C y (9) 

Y = (a y + iu)Y - \Y\ 2 Y + C x (10) 

where a x = a x — 2c x (1 + cos(7r/N)). The index j has been dropped in the Y array, since the 
array oscillates in-phase in steady-state. The subscript on C y (C x ) indicates that the term 
comes from the global coupling from the Y to the X array (from the X to the Y array). 
The constant components C x and C y in Eqs. © and (fTUj) are given by: 

C y = c xy Na l J 2 ; C x = c xy N[a x - 2c x (1 + cos(n /N))} 1 ' 2 . (11) 

Since c x < 0, C x > C y , it is clear that the X array has a stronger effect on the dynamics of 
the Y array than vice versa, even in the case of symmetric coupling: c xy = c yx . Equations of 
the form (0) and ()10|) have been studied previously in connection with forced Van der Pohl 
oscillators |2|. The main characteristic of these equations is that they execute a limit cycle 
below a certain critical value of the constants, C y , C x . Above the critical bifurcation value, 
the oscillations are damped out, so that the system reaches a stable equilibrium (a sink). 
There is also a narrow intermediate range (between the limit cycle and the sink region) which 
will not be discussed, but for a more detailed description see |2|. Equation JOJ) can be used to 
understand how ultraharmonic oscillations are created in the coupling term, c xy Ylk=i 1^1 
in Eq. (J3J). The constant component, C y , of the global coupling in Eq. (J0J), breaks the 
symmetry about the origin, making the amplitude of oscillation dependent on the phase of 
the oscillator. In phase-space, this corresponds to the shifting of the center of limit cycle 
from the origin. Now, all of the oscillators in the X array have amplitudes which oscillate 
at the frequency of the limit cycle, u. Since the phases of the oscillators in the X array 
are distrubuted in increments of 2tt/N, there are N peaks in amplitude over a single limit 
cycle oscillation. Adding all of the amplitudes together in the coupling term, c xy ^ k=1 \X k \ 
creates an ultraharmonic oscillation at frequency Nu. It follows that the global coupling 
term in the Y array will have a component oscillating at an ultraharmonic frequency, Nuo 
as well as a much larger constant component, C x . If this constant component, C x , is above 
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the bifurcation value, C x > Cb x , the Y array will be driven into the steady-state region of 
phase-space. It will then execute ultraharmonic oscillations with frequency Nu about that 
steady state. This can be better explained using averaging theory that applies when the 
system, given by Eq. (jTUJl is close to steady-state (when C x > Cb x ) and will be explored in 
the following section. 

The bifurcation values Cb x and Cby for C x and C y , respectively, occur when the the real 
part of the nullcline straightens out, leading to a single equilibrium, which is a steady state. 
Calculating the bifurcation value for Eqs. (0) and (JTUJ), 

where a = a y for Cb x and a = a x — 2c x (1 + cos(ir/N)) for Cb y - For a/uj > 1, Cb is a 
monotonically increasing function of a, so that the critical bifurcation for the X array is 
higher then for the Y array (for a x = a y ), since c x < 0. We are now in a position to give 
an explanation of the mechanism behind the onset of ultraharmonic oscillations in the Y 
array. As the symmetric global coupling constant, c xy = c yx , increases, a critical bifurcation 
value is reached in the constant component of the global coupling, whereby the Y array 
undergoes a bifurcation into a sink region of the phase-space, while the X array is still 
executing a limit cycle oscillation. After undergoing a bifurcation, the Y array oscillators 
are driven by the oscillatory ultraharmonic component of the coupling from the A-array 
and execute ultraharmonic oscillations about a steady state, which is determined by the 
constant component, C x . 

The constant component, C x , of the coupling term c xy X]fcLi l^fcl i n Eq. © is given by 

C x = c xy NA x (13) 

where A x is the average amplitude of the X array. If the amplitude of the X array is not 
substantially affected by the global coupling term, then A x can be approximated by Eq. (JHJ). 
Substituting Cb x from Eq. (JT2j) and A x from Eq. (jHJ) into Eq. (|T3*j) and solving for c xy , we 
get an approximation of the bifurcation value of the global coupling constant for the onset 
of stable ultraharmonic oscillations, 

b 2(a y /3f 2 + (3/16a y f 2 ^ 

° xy N[a x -2c x (l + cos(7r/N))Y/ 2 1 ] 

Figure (J2J) plots the bifurcation values of the global coupling constant, given by Eq. (fTlj) 
as a function of the diffusive coupling in the X array, c x . Note the excellent agreement 



between the analytically derived curve and the numerical simulation. The slight consistent 
under-estimate of c xy comes from using the unperturbed amplitude, A x , as given by Eq. 
(JSJ). As can be seen from the figure, however, the unperturbed amplitude, A x is an excellent 
approximation for the globally coupled system over a whole range of values. 

Figure (jHJ) plots the bifurcation values of the global coupling constant, given by Eq. (f"HJ) 
as a function of the parameter, a x . Higher values of a x , relative to a y lead to higher relative 
amplitude of oscillation of the X array, resulting in a lower bifurcation value of c b xy and a 
greater range of the global coupling constant, c b xy for generation of ultraharmonics. From 
Eq. (|14j) . for large values of the diffusive coupling, c x , the bifurcation value for the onset of 
ultraharmonic oscillations decreases as 1/ a/|c x |. Since Eq. (|T2|) is a monotonically increasing 
function of a, it is clear that the bifurcation value of the X array into a sink increases as a 
function of \c x \. It follows that as \c x \ increases, the range of the coupling constant c xy = c yx 
whereby stable ultraharmonics are induced in the F-array also increases. If the value of the 
global coupling is too high, then the X array undergoes a bifurcation, resulting in oscillation 
death. So, the coupling constant should be high enough, c xy > c xy so that the dynamics in 
Eq. bifurcate from limit cycle into a sink, but low enough so that the oscillators in the 
X array undergo a limit cycle oscillation at the natural frequency, uj. 



IV. AVERAGING THEORY AND THE CREATION OF ULTRAHARMONIC OS- 
CILLATIONS 

The mechanism behind the onset of ultraharmonics has been explained in the previous 
section. For small amplitude ultraharmonics, averaging theory can be used to prove that the 
frequency of these oscillations is the multiple of the limit-cycle frequency, Nu. Averaging is 
applicable to the systems of the form 

x = ef(x,t), xeH n (15) 

n 

where e is small and / is, T-periodic in t |2|. In our case, the periodic forcing comes from 
the ultraharmonic frequency, Nlj, in the periodic coupling term c xy Y2k=i l-^fcl- m system 
of the form given by Eq. ()15|). the relatively high frequency of periodic forcing constrasts 
with the slow evolution of the averaged system. Thus the averaging theory can be applied if 
the amplitude of ultraharmonic oscillations about the equilibrium given by solving Eq. (|T0j) 



S 



is small. As previously explained, the global coupling term can be broken up into a large 
constant component, C x , and the relatively small oscillatory component, C x . 



A' 



J2\Xk\ = C x + C x (t) (16) 



fe=i 



where 4 C C r Dividing Eq. (jHJ) by C x , we can now rewrite it as 

Y = F(Y)+C x (t)/C x (17) 



where variable Y has been rescaled, and F \ YJ is the right-hand side of Eq. (fTUj). scaled by 
C x . Near the equilibrium solution of Eq. (fTUj) . Eq. (JT7|) has the same form as Eq. (fT5|). since 

)oth C x (t)/C x and F are small. Applying the Averaging Theorem (see for example 

3]), we can thus solve for Y, 

Y = Y + W (18) 

where Y is the solution of the averaged system, given by Eqs. (jlUj) . and W is given by: 

dW 

Approximating C x to lowest order as a sinosoidal function at frequency Nu, we get an 
expression for Y as a function of C x , 

Y " ? + §t (20) 

Figure (jlj) compares numerically calculated amplitude of Y oscillation to the one given by 
Eq. (|20|). C x /Nuj. The agreement is better at lower ultraharmonic amplitudes, in accordance 
to the assumptions under which Eq. (J2Uj) was derived: C x {t)/C x <C 1. 

As can be seen from Eq. (I20J) . the amplitude of ultraharmonic oscillations is inversely 
proportional to the ultraharmonic frequency, Nu. We thus expect a degradation of ampli- 
tude as the number of oscillators in the array increases. The next section proposes a different 
coupling term that both achieves a higher amplitude of ultraharmonic oscillations, and does 
not suffer degradation in amplitude as iV increases. 

V. COUPLING FOR ACHIEVING HIGH FREQUENCY CONSTANT AMPLI- 
TUDE OSCILLATIONS 

Up to this point, we have focused on amplitude dependent global coupling, of the form 
given by Eqs. (j2J) and ©. For this form of coupling, the amplitude of ultra-harmonic oscilla- 
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tions is low and falls as N increases. The amplitude of ultraharmonics cannot be significantly 
increased by increasing the mutual coupling strength, since high mutual coupling leads to 
a bifurcation of the X array and oscillator death. From the previous sections, it should be 
clear that any form of coupling, which is a periodic function of some frequency, uif, induces 
oscillations at that frequency (for sufficiently low amplitude) as long as the constant part 
of the coupling, C x , is above the bifurcation value: C x > Cb x - It should therefore be pos- 
sible to better control the amplitude of ultra-harmonic oscillations and even the frequency 
by the choice of the coupling function. In this section, a form of coupling is derived that 
induces ultraharmonic oscillations in the Y array at frequency 2uN (rather then uN), with 
an amplitude that stays relatively constant as iV increases. N used will be odd, since as 
previously discussed, only N odd arrays have phases in increments of 2ir/N in the presence 
of diffusive coupling. 

As shown in Section II, in the absence of inter-array coupling, c yx = 0, the X array 
oscillates as a collection of simple harmonic oscillators, at an amplitude given by Eq. (JHJ), 
and with a phase difference of ir + n/N between nearest neighbors and 2-7T /N jumping over 
a neighbor. The function describing each oscillator in a steady-state is given by: 

Xj r = A x COS {ojt — (fij) (21) 

where Xj r denotes the real part of the jth oscillator, Xj, and A x is the amplitude given by 
Eq. (|HJ). The phase, <fij is given by 

& = (;-l)(7r + 7r/iV) (22) 

where Eq. (jSJ) was used. Since {Xj r } is a collection of sinosoidal functions with phases 
distributed in equal increments over the interval [0, 2tt], an ultra-harmonic coupling function, 
c xy Ylk=j can thus be created by simply summing the real part of the oscillation, 

c xy J2k=i — J2f=i l-^jrl- I n this case however, the constant component of the coupling, 

C x increases substantially as iV increases, and the relative amplitude of the ultra-harmonic 
component, C x drops substantially with N (where Ylj=i \^jr\ = C x + C x ). Figure El shows 
this effect for N = 3 and N = 5. It is possible to increase the amplitude of the oscillatory 
component of the coupling and decrease C x by taking the sum of some power of Xj r , 

N N 

c xy Y,9(Xj) = c xy Nj2\X jr \ n (23) 

3=1 3=1 
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where n is an integer yet to be determined, and Xj r has been normalized (divided by A x ): 

Xj r = COS {ut — (24) 

where <frj is given by Eq. (|22|). The factor of N multiplying Eq. ()23)1 is there to compensate 
for the fall of ultraharmonic amplitude as the frequency increases (see Section IV). We need 
to find an expression for n as a function of iV and C x , n(N, C x ), such that the amplitude of 
oscillation stays fairly constant as iV increases. 

From Eq. (|24jh it is clear that increasing n (for even n) in Eq. will lead to sharper, 
more narrow peaks centered around cut — {<f)j} = 0, 7r, with j running between 1 and N. 
Thus with increasing n, there will be less overlap between neighboring peaks in Eq. (|24j) . 
leading to a lower value of C x , the constant component of the coupling, (see Figure IBJ). In 
the limit as n — > oo, the sum in Eq. becomes, as a function of t, 




where S is the delta function. In steady state, and for n-even, the normalized coupling from 
the X array as n — > oo is given by a series of equally-spaced spikes, occuring at frequency 
2nu and of amplitude c xy . Figure |H1 shows Ylk=i \Xkr\ n f° r n — 60. We can already see 
sharply defined spikes that approach a sum of spikes of a unit amplitude as n — > oo. Eq. 
(f23j) has frequency of 2Nu rather then Nu, the frequency generated when the amplitudes 
are added. This happens because each |Xfc r | n has 2 spikes over the time interval [0, 2n/u], so 
that summing over N oscillators leads to a waveform of frequency 2uN, rather then uN, as 
for amplitude coupled arrays (see Figured). Since neighboring peaks are separated by n/N, 
the two nearest spikes intersect at a phase difference of tc/2N from the top of each peak. 
Assuming that the exponent, n, is sufficiently large so that only nearest neighboring spikes 
have significant overlap (see Fig. EJ), we are led to the following equation for the constant 
component of the coupling from the X-array to the V-array: 

C x = 2N-c xy [cos(ir/2N)] n , (26) 

where the exponent, n is a function of N, such that cos (tt /2N) stays constant as N increases. 
Therefore C x increases linearly with N, a situation similar to the amplitude coupled arrays, 
where the constant component of the coupling also increases linearly as N, for large N 
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(see Eq. ([T3])). The amplitude of the oscillatory component of the coupling is found by 
subtracting the constant component, C x , given by Eq. (|2l)j) from the peak of the normalized 
oscillation, which occurs at one, 

max(C x ) = c xy N (1 - 2[cos (n/2N)} n ) . (27) 

Thus, unlike amplitude coupled arrays, the oscillatory component of the coupling also in- 
creases linearly with N, which prevents the degradation of amplitude seen in the amplitude 
coupled case. Using Eq. (J2*U|) with 2Nuj in the denominator (since that is the frequency of 
the drive), the amplitude of ultraharmonic oscillations can be approximated as 

A v w^(l-2[co S (7r/2iV)r) (28) 

and should therefore stay fairly constant as N increases if [cos (ir/2N)] n ~ const for large 
N. For large N, we can approximate the cosine term as 

2N\ 

where the first approximation came from taking the first two terms of a series expansion of 
a cosine function and the second approximation came from expanding the exponential in a 
series, taking the first two terms and comparing them to 1 — (7r/2N) 2 . Thus in order for 
A y in Eq. (|28|) to stay fairly constant as N increases, for N large, we need n oc N 2 . When 
taking the power, n should be rounded to the nearest even number. Choosing 

n = KN 2 - KN 2 mod2 if KN 2 mod2 < 1 (30) 



cos{ix/2N) w 1 - (n/2N) 2 w exp [ — ^— ) (29) 



n = KN 2 - KN 2 mod2 + 2 if KN 2 mod2 > 1. (31) 

The above equations ensure that n is rounded to the nearest even number. K is some 
constant, which can be chosen to achieve a desired value of C x (for example K = 2/3 is a 
good choice). Using Eqs. (f28|) - (j3~Tf we obtain an approximate expression for the amplitude 
of ultraharmonic oscillations when iV is large: 

From the above equation, the coupling strength, c xy and K can be chosen to obtain the de- 
sired amplitude of ultraharmonics, with a maximum possible amplitude being A y w c xy /2ui. 
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Figure (JHJ) shows the amplitude of ultraharmonics, A y , as a function of N for both the 
drive coupling given in Eq. ([23)1 . with n given by Eqs. ()30|) and (J3*TJ) and mutual coupling 
analyzed in the previous section. The parameters were chosen such that the constant com- 
ponent of the coupling, C x is the same in both cases. The solid line for the drive dependent 
coupling is the result of connecting the data points for the odd N values between N = 3 
and AT = 31. The amplitude of ultraharmonics asymptotes to a constant value of around 
0.033, close to the value predicted by using Eq. (J3*2*J) of around 0.028. It is clear that the 
amplitude of ultraharmonics generated by the drive is significatly higher than that by the 
mutual coupling, given by Eqs. (J2J) and (JHJ). In addition, the drive coupling does not suffer 
from oscillator death, thus the amplitude will stay constant while N increases. In contrast, 
increasing N for mutual coupling leads to oscillator death, unless the coupling constants 
c X y,c yx are also decreased (leading to an even further fall in amplitude). In Figure (JHJ), 
oscillator death occurs at just N = 9, for the coupling used. Thus the drive has many 
advantages, such as a relatively high, easily controllable amplitude that does not degrade 
with an increase in N and does not suffer from oscillator death, which occurs in mutually 
coupled arrays. However, since n increases as N 2 , high ultraharmonics, Nu, require rather 
high powers in the coupling function. 



VI. CONCLUSION 



The mechanism behind the generation of ultra-harmonic oscillations in two mutually 
coupled arrays of limit cycle oscillators was analyzed. These ultra-harmonic oscillations were 
shown to occur as a result of a bifurcation that results when the coupling from the X to the Y 
array exceeds a certain value. This coupling consists of a large constant component, since the 
coupling is amplitude-dependent, and a smaller ultraharmonic oscillatory component that 
results from the breaking of symmetry in the X-array due to coupling from the K-array. This 
smaller oscillatory component of the coupling induces the Y array to oscillate at an ultra- 
harmonic frequency around the equilibrium determined by the large constant component of 
the coupling. It was also shown that in the case of amplitude-dependent coupling, the ultra- 
harmonic oscillation is the result of a mutual interaction between the two arrays (rather 
than a master-slave system), since the coupling from the Y to the X array is necessary to 
induce an oscillation in the otherwise conserved amplitude of X oscillators. 
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The range of inter-array coupling constants whereby ultra-harmonic oscillations are pro- 
duced was derived for symmetric c xy = c yx coupling, but can be generalized to non-symmetric 
coupling, where c xy ^ c yx . For symmetric coupling, this range depends on the limit cycle 
frequency, the amplitude of oscillation and the strength of nearest-neighbor coupling within 
the X array. The allowable range of inter-array coupling constant increases with an increase 
in strength of the diffusive coupling, c x , in the out-of-phase coupled X array Thus higher 
absolute values of c x lead to a greater possible range of values of the inter-array coupling 
constant whereby ultra-harmonic oscillations are created. The Y array diffusive coupling 
strength, c y , on the other hand, does not affect the conditions for creation of stable ultra- 
harmonics, since the diffusive coupling term drops out in steady-state for in-phase coupled 
arrays. The derived bifurcation values for the inter-array coupling constant agree well with 
numerical simulation, and can be used to tune the value of the coupling constant to control 
the amplitude of the ultra-harmonic oscillations. 

For achieving better-controlled, higher amplitude ultra-harmonic oscillations, another 
form of coupling was suggested. This one-way coupling has the advantage of achieving 
higher amplitude ultra-harmonic oscillations that do not fall in amplitude as the number 
of oscillators, and therefore the ultra-harmonic frequency, increases. It also does not suffer 
from oscillator death, which puts an upper limit on the strength of symmetric coupling that 
can be used in amplitude coupled arrays. The suggested form of coupling, however, requires 
increasingly more complicated forms of the coupling function as iV increases and may be 
more difficult to implement experimentally. 

Finally, though the methods of analysis here were applied to oscillator of Stuart-Landau 
type, they may be applied to various applications of interest, where both frequency and 
power control are required. Such examples of the stabilization of in-phase arrays occur in 
such areas as electronic circuits for radar [3| , phase locked nano-scale magnets used for 
microwave sources arrays used for terahertz 
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{Xj}: Out-of-phase 
diffusive coupling 




{Yj}: In-phase 
diffusive coupling 



global all- to-all 
inter-array coupling 



N=5 




FIG. 1: Two diffusively coupled arrays coupled to each other via global coupling, for N = 5. The 
X array has out-of-phase diffusive coupling, and the Y array has in-phase diffusive coupling. The 
two arrays are globally coupled to each other. 
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FIG. 2: Bifurcation values of the global coupling constant, c xy as a function of out-of-phase diffusive 
coupling, c x . The solid line is a plot of the analytically derived Eq. (|14[) . The numerical values are 
also shown on the graph 
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FIG. 3: Bifurcation values of the global coupling constant, c xy as a function of amplitude, a. 
Both the analytically derived dependence and the numerics are shown 
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FIG. 5: Top: Y2j=i \Xjr\ f° r N = 3. Bottom: J2j=i \^jr\ f° r N = 5. In both cases a = 1, 
u) = 1/2, = 0. As increases, the constant component for this form of coupling increases, 
while the relative value of the oscillatory component falls. 
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FIG. 6: Drive given by Eq. (j2.3|) for high n (n = 60). The peaks narrow substantially, with less 
overlap between neighboring peaks as n increases. 
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FIG. 7: Oscillations with coupling given by Eq. (|23[). Top: the coupling term, c xy N Ylf=i l-^-ji 
with C x = 0.8 and C x = 0.8. Bottom: The ultraharmonic oscillations induced in the Y array 
the coupling shown at the top. a = 1, u) = 1/2, c xy = 0.36, c yx = 0, N = 3. 
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FIG. 8: Amplitude of ultraharmonic oscillation, A y , with the drive coupling, as a function of N, 
K = 2/3. The lower three points, marked by crosses, correspond to amplitude of ultraharmonics 
with the standard amplitude coupling. The constant part of the coupling, C x is the same for both 
cases, c xy = c yx = .1134, a x = a y = 1, u = 1/2 
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